Determining optimal cloud cover
We can now analyze these search results so that we only download the tiles we need. One paramount parameter is the amount of clouds identified in the imagery. To define the threshold we defined min_SAFE as the minimum number of sensed tiles across all tiles for a given cloud threshold and sum_SAFE as the total of sensed tiles for a given cloud threshold.
search_stats <- search_results %>% dplyr::group_by(id_tile) %>% dplyr::summarise(dplyr::across(dplyr::contains("clouds"), list(min = min, mean = mean))) %>% dplyr::summarise(clouds_min = min(clouds_min), clouds_mean = mean(clouds_mean))
clouds_stats <- data.frame(
clouds_threshold = seq(ceiling(search_stats$clouds_min), ceiling(search_stats$clouds_mean), by = 0.1)
) %>%
dplyr::mutate(
min_SAFE = sapply(clouds_threshold, function(clouds_threshold){
search_results %>%
dplyr::group_by(id_tile) %>%
dplyr::group_map(~ {dplyr::filter(., clouds <= clouds_threshold) %>% nrow()}) %>%
unlist() %>%
min()
}),
sum_SAFE = sapply(clouds_threshold, function(clouds_threshold){
search_results %>%
dplyr::group_by(id_tile) %>%
dplyr::group_map(~ {dplyr::filter(., clouds <= clouds_threshold) %>% nrow()}) %>%
unlist() %>%
sum()
})
)The graph below shows the trade-off between incorporating more cloud: the minimum number of sensed tiles goes up, but so does the total number of sensed tiles. Each tile is roughly 1 GB.
p <- ggplot2::ggplot(clouds_stats %>% reshape2::melt(id.vars = c("clouds_threshold")), ggplot2::aes(x = clouds_threshold, y = value)) +
ggplot2::geom_point() +
ggplot2::facet_wrap(. ~ variable, scales = "free") +
ggpubr::theme_pubclean()
plotly::ggplotly(p)From these two graphs, we select the following threshold:
sensing_limit <- 30
clouds_stats %>% dplyr::filter(min_SAFE > sensing_limit) %>% head(1)
#> clouds_threshold min_SAFE sum_SAFE
#> 1 5.1 31 4166
clouds_threshold <- clouds_stats %>% dplyr::filter(min_SAFE > sensing_limit) %>% head(1) %>% dplyr::pull(clouds_threshold)which is appropriate as the cloud percentage remains low, at least 30 tiles are sensed for any tile, and less than 3 TB of SAFE data have to be downloaded.
The selection of the product for downloading is then done by:
selected_search_results <- search_results %>% dplyr::filter(clouds <= clouds_threshold)
dim(selected_search_results)
#> [1] 4166 12Sensing dates
Additionally, we can visualize the sensing dates with the following codes:
p <- ggplot2::ggplot(selected_search_results, ggplot2::aes(x = sensing_datetime, y = clouds, color = id_tile)) +
ggplot2::geom_point() +
ggpubr::theme_pubclean()
plotly::ggplotly(p)Visual inspection of the above graph tile by tile does not indicate that there is a significant skew towards a given period for any tile. Example of a tile with discontinuous record is 10TGK. The tile with the lowest average cloud cover is 10SEF.
Number of sensed tiles
The following graph shows the distribution of the number of sensed tiles:
count_df <- data.frame(count = selected_search_results %>%
dplyr::group_by(id_tile) %>%
dplyr::group_map(~ nrow(.)) %>%
unlist()
)
ggplot2::ggplot(count_df, ggplot2::aes(x = count)) +
ggplot2::geom_density() +
ggplot2::labs(x = "number of sensed tiles") +
ggpubr::theme_pubclean()